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Identification of the system parameters of a randomly excited 
structure may be treated using a variety of statistical tech- 
niques. Of all these techniques, the Random Decrement is uni- 
que in that it provides the homogeneous component of the 
system response. Using this quality, a system identification 
technique was developed based on a least— squares fit of the 
signatures to estimate the mass, damping, and stiffness 
matrices of a linear randomly excited system. In this part of 
the paper the mathematics of the technique is presented in 
addition to the results of computer simulations conducted to 
demonstrate the prediction of the response of the system and 
the random forcing function initially introduced to excite 
the system. Part II of the paper presents the results of an 
experiment conducted on an offshore platform scale model to 
verify the validity of the technique and to demonstrate its 
application in damage detection. 


INTRODUCTION 

In general, all system identification techniques begin by assuming a form for 
the equations describing the system, then attempt to identify the unknown parameters 
in that assumed system through prior knowledge of the actual response, and sometimes 
the input as well. For linear systems, the identification process could be conducted 
in two different ways depending on the available information. To describe the two 
methods, consider a multidegree-of-f reedom system having the following set of dif- 
ferential equations: 


[M] X + [K] X = F 


( 1 ) 


where X is the response vector and F is the input loading vector. The response vec- 
tor may further be viewed as the sura of the homogeneous solution vector X. and the 
particular solution vector X , i.e. “ 

P 


X 


*h + X p 


( 2 ) 
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The first way of identifying the system in Equation (1) is to measure the input 
into the system, F, and the response of the system, X. Then through the use of a 
curve fitting technique, matrices [M] and [K] may be identified [14,15]. The disad- 
vantage of this method is that the input must either be of a type that could be 
measured, or it must be deliberately introduced into the system. Furthermore, if the 
system has N degrees-of-f reedom and only M locations are monitored, where H < N, 
these monitored locations must be selected specifically to include all external 
loads into the system. Therefore, this method is impractical in applications where 
the system is naturally excited, such as offshore structures impacted by wave 
motion, and flight vehicles excited by turbulent air flow. 

The other approach for identifying the system parameters is through the use of 
the system response only [2,7,8,9,12]. This leads to the identification of the 
eigenvalues and eigenvectors of Equation (1) as opposed to the mass and stiffness 
matrices. 

Substituting Equation (2) into Equation (1) yields 


[M] + x p ) + [K] ^ + x p ) - r 

which may be separated into two independent equations, namely 
[M] X p + [K] X p = F 
[M] + [K] Xb = {0} 


Equations (3) and (4) indicate that if F is not known, the system parameters may 
not be identified unless the homogeneous and particular solutions are separated. In 
practice, however, if the input spectral density is relatively flat over the range 
of the system frequencies, the ratios of the responses x^(t) at the different loca 
tions are taken in the frequency domain thus yielding the eigenvalues and the eigen- 
vectors of Equation (4). If, on the other hand, the input spectral density has some 
mild fluctuations over the frequency range of the system, then taking the ratios of 
the responses could yield erroneous eigenvectors. In addition, if the system exhi 
bits some damping, the modal damping ratios may not be identified correctly. This is 
mainly due to the fact that the frequency content of the particular component of the 
response may vary considerably at different locations in the system. 

These problems may be overcome by employing the Random Decrement (Randomdec) 
and cross-Random Decrement (cross-Randomdec) techniques [3,4,5,6,10,11,13]. 

Given the response vector X, the Randomdec and corresponding cross-Randomdec signa- 
tures are calculated. Based on the results obtained in reference [1], the signatures 
should be interpreted as the homogeneous components of the response, namely X^. 

With this being the case, the eigenvalues and eigenvectors of Equation (4) may be 
estimated accurately. The biggest advantage of this method is that the response need 
not be measured at the points where external loads are applied. Furthermore, modal 
damping may be evaluated accurately from the Randomdec signatures either in the time 
domain using the logarithmic decrement or in the frequency domain using the half 
power point (curve fitting could generally be used for close modes in both domains). 

Although modal parameters might be sufficient for many applications, it is some- 
times desirable, if not necesssary to have information reflecting the actual mass, 
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damping, and stiffness of the system. This could be the case in systems where 
damping is not proportional, or when the model is required for damage detection in 
which elements in the original matrices pertain to actual locations in the system. A 
technique is therefore proposed to estimate the [M], [C], [K] matrices of a linear 
system with the use of the Random Decrement technique. 


PROPOSED SYSTEM IDENTIFICATION TECHNIQUE 
Consider the linear set of equations 


[M] X + [C] X + [K] X = F 


(5) 


where [M] and [K] 
symmetric damping 

are real 
matrix. 

symmetric matrics 
Introducing matrix 

and [C] 

H . and 
Pij 

is a nonproportional 
vector Z pj , where 
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Equations (5) may be rewritten in the form 


M 


l 

j=l 


l H 

p=l 


pi j 


Z . 
PJ 


i = 1,2, 


,M 


( 6 ) 


where M is the number of degrees - of "freedom in the system, and F. symbolizes ele- 
ment i of vector F. 1 

If Z and F are composed of N discrete points in time, there should exist one 
equation 15 similar to Equation (6) for every point k, where k = 1,2,...,N . Therefore, 
for one time step k 


M 3 

V V H Z = F 

p =! Pij PJk ik 


i = 1,2, ... ,M 
k = 1,2 N 


( 7 ) 


To identify the three matrices H , a least squares scheme will be employed to 
obtain the best estimate for H through minimizing the difference between the left 
and right sides of Equation (7). Therefore, defining an error index e., for each 
equation i at every time step k. Equation (7) may be rearranged as follows : 
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where H . . is the best estimate for H . .. 

Add?n| the sum of the squares of Equation (8), the total error E in the system 

may be defined as 
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(9) 


To 

to all 


minimize the error with respect to H . , the slope 
the unknown parameters must be set tB ^ero, i.e. 


of Equation (9) relative 


9E 


9H 


pi j 


0 


P = 1,2,3 
i = 1,2, ... ,M 
j = 1,2, ... ,M 


Therefore, taking the partial derivative of E with respect to H ^ 

and noting that ft . = H (symmetry condition), the following two equations result 
pi J PJi 


9E 


9H 


pi j 


- 2 JLv [ F * - 1 X ^ v 1 ' ° 


(10a) 


for i = j , and 


9E 


9H 


pi j 


2 J^ Z pjk [ F ik 'J x H qin Z qnk ] 


+ Z pik I F jk ~X X \in Z qnk H = ° 


(10b) 


for i * j . 

Equations (10a) and (10b) form a set of M(M + 1)3/2 linear simultaneous 
equations which may be solved either in closed form or by iteration. To use the 
latter, the partial derivatives must be taken of Equations (10a) and (10b) with 
respect to all the unknown parameters resulting in a constant Jacobian matrix, 
name ly 
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for i = j , and 
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N 


\l 1 V V + V V i 


for i * j . 

To obtain a solution in closed form, Equations (10a) and (10b) should be rewrit- 
ten in the form 


N M 3 

tJl n-1 q=l Hqin Z P ik Zqnk 


N 

y z f 

k ; L p^ ik 


(Ha) 


for i = j , and 
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J pjk ik 


Z F 1 
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(Hb) 


i° r i * J> where the right side of Equations (11a) and (lib) contains the constant 
terms while the left side is a constant coefficient matrix premultiplied by a vector 
composed of the M(M + 1)3/2 unknown system parameters. This may be represented as 
follows: 


EG] h = F 


( 12 ) 


These equations can now be solved using any conventional linear equation solver. 
It is apparent by inspection^ Equation (12) that when the Random 

Decrement technique is applied, F reduces to zero , thus resulting in 


[G] h = {0} 

or, in expanded form 
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(13) 
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where g = M(M + 1)3/2. Equations (13) are in homogeneous form and therefore do not 
possess a unique solution for h* 

On the other hand, if the rank of [G] is g - 1, and if one of the unknown para- 
meters in h were indeed known, then a unique solution would exist. This may be 
proved by rearranging Equation (13) as follows : 
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(14) 


where h was assumed to be the known parameter for convenience. Equation (14) may be 
set up for any known h by simply eliminating the i th row of [G] (since the partial 
derivative with respect to a constant is zero) then moving the i th column to the 
right side and multiplying it by the scalar -h . The only condition for selecting 
h^ is that the i^h row must increase the rank of the system. With the aforemen- 
tioned steps adhered to. Equation (14) will definitely have a unique solution. 
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If the rank of matrix [G] is less than 3 - 1 » the solution of Equation (13) 
becomes more difficult to obtain. The linear dependency of the equations forming 
matrix [G] will require a greater effort to be eliminated, and, more parameters of 
vector h will have to be assumed. Therefore, if the rank is 3 “ 2, two equations 
from Equations (13) must be removed and two parameters of vector h must be assumed. 
Calculating two or more parameters in many physical systems is extremely difficult, 
if not impossible. The examples used in the remainder of this paper are all of 
systems with a rank of 3 “ 1* 

The question to be raised then is how can the value of one of the parameters be 
known ? 

Since h is composed of elements from the mass, damping, and stiff-ness matrices, 
it is physically immeasurable. Therefore, the best alternative would be to obtain a 
good estimate. The requirement imposed on the accuracy of the estimate depends 

greatly on the application for which the model is needed. If the model is to be used 

for damage detection, then knowing the absolute values of the [M], [C], [K] matrices 

is not essential, but rather, the ratio of the values of the elements at different 
instances in time. Therefore, assuming the order of magnitude of one of the elements 
in h should suffice as long as the same value is used every time. On the other hand, 
if the model is needed for simulation purposes, where one wishes to study the effect 
of different loading conditions on the response of the system, the absolute values 
are needed. In this case, the estimate for h may be obtained from either a finite 
element model, or by performing a simple static test on the system (if physically 
possible). The latter is carried out by applying a known static load at a point in 

the system and measuring the deflection at the same point or at any other point. 

These two values may then be used to scale the entire h vector (h must be already 
calculated by assuming the order of magnitude of h^). The scaling procedure may be 
described by noting the static component of Equation (5), namely 

[K] X = F (15) 


If [K] is separated into a scalar e multiplying a matrix [K], 

where [K] is the estimated stiffness matrix using the proposed method, 
then Equation (15) becomes 

e [K] X = F (16) 

Since F and [K] are known (F is the force applied during the static test) the 
system of equations 


[K] X - F 


may be solved to give X , where 


X = e X 
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But, the response X at point i in the system was measured, therefore 
e may be obtained by simply taking the ratio X^/X^ . The value of e 
may then be multiplied by the estimated mass, damping, and stiffness 
matrices [M] , [C] , and [K] respectively to obtain [M], [C], and [K]. 


TESTS OF PROPOSED TECHNIQUE BY COMPUTER SIMULATION 

To demonstrate the procedure for applying the proposed system identification 
technique, and to test its accuracy, a six degree-of-f reedom model of a cantilever 
beam is used. Figure 1 shows the beam with the locations of the six points where the 
response is monitored. Using a finite element model (Reference [15]), the system 
matrices in Equation (5) were obtained. 
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Figure 1 - Locations of monitored points on the cantilever beam 


To distinguish between errors introduced by the system identification technique 
and the errors in the Randomdec signatures, two different scenarios were conducted. 
The technique was initially tested using exact free-decay response curves of the 
system and then tested using actual Randomdec signatures obtained from the random 
response of the system. 

Case I: Exact Free-Decay Curves 

A set of initial conditions was arbitrarily chosen for the six locations on the 

beam. Equation (5) was then solved numerically and the response vectors X, X, and X 
were recorded. The time step size was selected to insure that at least seven points 
were needed to construct one cycle of the highest frequency in the system. Solving 
for the eigen-values of Equation (5) the undamped natural frequencies of the system 
were found to be 

f, = 25.908 Hz f = 886.770 Hz 

1 4 

f 2 = 159.251 Hz f 5 = 1445.039 Hz 

f. = 445.307 Hz f = 2114.293 Hz 

J o 
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thus resulting in a step size of AT = 0*00007 seconds. Furthermore, the number of 
time steps had to be selected to cover at least one cycle of the lowest frequency in 
the system. Therefore, 600 points were used* 

One of the tests for the technique is its repeatability relative to different 
fixed parameters, i*e. as elements h^ through h^ are in turn fixed. Since diagonal 
terms are more reliable than off-diagonal terms in the model, this test was strictly 
confined to the diagonal elements. 

Equation (14) was solved eighteen times, each time fixing one of the diagonal 
elements in the [M], [C], [K] matrices. The fixed value was always taken as the 
actual value of the element to avoid scaling the h vector. After every evaluation 
of the h vector, the errors occurring in the diagonal terms, relative to their 
actual values, were calculated. Calculations were made of the errors in the mass, 
damping, and stiffness matrices as each diagonal element in the mass matrix was 
fixed. Similar results were obtained as each diagonal element in the damping and 
stiffness matrices was fixed, respectively. In all the cases, the average error 
occurring at every point on the beam was evaluated. 

Following a careful inspection of the errors, two conclusions were made, viz. 
fixing the stiffness matrix gives better estimates of vector h than fixing the mass 
and damping matrices, and , due to the large variation in the error resulting from 
fixing different elements in the same matrix, the avarage value is probably a more 
consistent estimate. Therefore, based on these conclusions, and noting that fixing 
the stiffness matrix elements is the only case in which the average gives a better 
overall estimate than the individual estimates, it is further concluded that the 
best approach for estimating vector h is to take the average of the h vectors 
obtained by individually fixing the stiffness diagonal elements. In doing so, the h 
vectors must be scaled independently before the averaging process is carried out. 
Using this procedure, the system matrices were identified. These matrices were not 
identical to the actual system matrices; their validity to represent the system was 
checked by comparing the reponse to the same input. Therefore, a random input vec- 
tor F with constant spectral density was simulated on the computer and used as input 
into Equation (5). This equation was solved using the actual and the estimated 
system matrices. Results were obtained for all the response points. Comparison of 
the actual and predicted responses at point 1 is shown in Figure 2. The responses 
of the two systems compare favourably. 



Figure 2 - Comparison of actual and predicted responses at point 1 
(system identification from actual free-decay response) 
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Case II: Free-Decay Curves From Randomdec Signatures 

To further evaluate the overall accuracy of the system identification technique, 
the technique had to be tested by incorporating Randomdec and cross-Randomdec signa- 
tures as opposed to actual free-decay response curves. 

Two purely uncorrelated stationary, Gaussian, random records were used as input 
forces at points 2 and 4 on the cantilever beam. Equation (5) was solved and the 
response vector X recorded. Due to its high frequency content, station 6 was used as 
the triggering station with a trigger level of 0.0075. Five cross-Randomdec signa- 
tures (at stations 1 to 5), and one Randomdec signature (at station 6) were obtained 
for 600 lag points and 500 averaged segments. The cross-Randomdec signature for 
station 1 is shown in Figure 3. The first and second derivatives were then calcu- 
lated for t^e six signatures using a finite difference scheme with an error on the 
order of AT . 



Figure 3 - Cross-Randomdec signature of time record at point 1 


The procedure recommended in the previous section was used. Equation (14) was 
solved six times, each time fixing one of the diagonal elements in the stiffness 
matrix. was fixed at 1000, was fixed at 5000, and through K we re 

fixed at 10000. In each case, a load of 10 was applied at point 1 and, using the 
corresponding estimated stiffness matrix, the deflection at point 1 was calculated 
(since point 1 is the free end of the cantilever beam). The same procedure was 
followed using the actual stiffness matrix. The six h vectors were then scaled 
following the procedure outlined in the previous section, averaged, and rearranged 
in matrix form. 

Again these matrices were not identical to the actual system matrices. To test 
their validity as a simulation tool, the estimated matrices were substituted into 
Equation (5). The force vector used to obtain the random records, from which the 
signatures were evaluated, was used as the input. Comparison of the calculated 
response versus the actual system response was obtained for all the stations. 
Results at stations 1 and 2 are shown in Figures 4 and 5. Once more, the results 
compare favourably . 
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Figure 4 - Comparison of actual and predicted responses at point 1 
(system identif icatin from Randomdec signatures) 



TIME (SEC) 

Figure 5 - Comparison of actual and predicted responses at point 2 
(system identification from Randomdec signatures) 


ESTIMATING THE INPUT INTO THE SYSTEM 

An interesting application arises from the system identification technique by 
observing Equation (5). After estimating the mass, damping, and stiffness matrices, 
if the response vector X and its derivatives are substituted back into Equation (5), 
the outcome should be the force vector F. 

Therefore, in a real application, the random response would be measured at 
several locations in the system. The signatures would then be obtained, their deri- 
vatives calculated, and Equation (14) used to estimate the [M], [C], [K] matrices 
following the procedure outlined in the previuos section. A direct substitution of 
the measured response and its derivatives into Equation (5) with the estimated 
matrices would result in a vector similar to the input vector. If the estimated 
mass, damping, and stiffness matrices were not scaled, the outcome of Equation (5) 
should be a scaled version of the input vector. 
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To demonstrate this approach, the example provided in the previous section was 
used. Results were obtained of the estimated input vector F versus the actual input 
Since the forces were originally applied at locations 2 and 4 on the beam, the for- 
ces at locations 1, 3, 5, and 6 should be zero. These Results indicate that at the 
points where the loads were applied, the estimated input functions formed good 
approximations (see Figures 6 and 7). As for the unloaded points, the technique pre 
dieted forcing functions with relatively small magnitudes in comparison to the 
loaded points (see Figure 8). 



Figure 6 - Comparison of actual and predicted force records 
(record for location' 2 - unfiltered response) 



Figure 7 - Comparison of actual and predicted force records 
(record for location 4 - unfiltered response) 
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Figure 8 - Comparison of actual and predicted force records 
(records for locations 1,3,5, and 6 - unfiltered) 


Care must be excercised when employing filters in this technique. If the signa- 
tures are obtained after filtering the response record, the response vector X and 
its derivatives must also be fitered before substitution into Equation (5) . If this 
procedure is not followed, the estimated force vector F will also include the 
filtered modes of the system. This may be demonstrated using the same example. 
Studying the Fourier magnitude spectrums of the responses at points 1 and 6 and 
their derivatives , it was apparent that the lowest mode is quite dominant. 
Therefore, employing a high-pass filter at 80 Hz for the signatures, but not the 
response vector X, the resultant estimated inputs at points 2 and 4 are shown in 
Figures 9 and 10. The low mode is quite apparent in the estimated input records. 

This problem could also occur when calculating signatures from velocity and acce- 
leration records since they usually tend to include a larger density of the higher 
frequency inodes. 

In addition to the aforementioned effects, the frequency content of the response 
records is a major cause for the dissimilarity between the identified system matri— 
ces and the actual matrices. This may be explained by considering the six degrees— 
of^freedom cantilever beam. If the beam were excited by a random force with a 
band-limited frequency range, where, for the sake of example, this range included 
only the lowest three modes, then the identified matrices will possess information 
concerning these three modes only. This would mean that the identified [M], [C], 
and [K] matrices would definitely be different from the actual matrices. Therefore, 
subjecting the predicted system to the actual input vector should result in a 
response very similar to the measured system response, whereas a force rich in the 
higher frequencies would yield different results. Hence, using the proposed tech- 
nique for damage detection requires the frequency content of the input force to 
always be the same. 
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CONCLUSIONS 


A system identification technique was proposed based on a leastsquares fit of 
Randomdec and cross-Randomdec signatures to identifiy the mass, damping, and stiff- 
ness matrices of a linear multidegree-of-f reedom system. Computer simulations 
carried out for a discretized finite element model of a cantilever beam proved the 
technique to be quite effective in predicting the response of the beam for a given 
frequency range of excitation. Furthermore, the proposed technique was demonstrated 
to be successful in predicting the random forcing function initially introduced to 
excite the system. The results of the simulation clearly indicated the importance of 
filtering the response of the beam and the effect it may have on the identified 
system. 
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NOMENCLATURE 
flexibility matrix 

damping matrix of multiple D.O.F. system 
unsealed identified damping matrix 
element ij of damping matrix [C] 
input loading vector 
forcing function applied at point i 
system matrix to be identified 
estimated system matrix 

vector containing system parameters to be identified 
element i of vector h 

stiffness matrix of multiple D.O.F. system 
unsealed identified stiffness matrix 
element ij of stiffness matrix [K] 
mass matrix of multiple D.O.F system 
unsealed identified mass matrix 
element ij of mass matrix [M] 

matrix containing ratio of flexibility matrices 
element ij of matrix [R] 
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t 

X 


*p 

x i 


time variable 

response vector of multiple D.O.F. system 
homogeneous response vector 
particular response vector 
system response at point i 


x(t),x(t),x(t) 


e 

e 


position, velocity, and acceleration of variable x 
measured system response 

number of system parameters to be identified 
scaling factor of identified system 
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